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Abstract.  The  general  form  of  Maxwell’s  velocity  slip  boundary  condition  for  rarefied  gas  flows  depends  on  both  the 
geometry  of  the  surface  and  the  constitutive  relations  used  to  relate  the  viscous  stress  to  rate  of  strain.  The  dependence  on 
geometry  is  often  overlooked  in  current  rarefied  flow  ealculations,  and  the  generality  of  the  constitutive  dependence 
means  the  eondition  can  also  be  usefully  applied  in  regions  where  the  Navier-Stokes  equations  fail,  e.g.  rarefied  flows 
close  to  surfaces.  In  this  paper  we  give  examples  illustrating  the  importance  of  both  these  dependencies  and  show, 
therefore,  that  implementing  the  general  Maxwell  eondition  produces  substantially  different  results  to  conventional 
implementations  of  the  condition.  Finally,  we  also  investigate  a  common  numerical  instability  associated  with  Maxwell’s 
boundary  condition,  and  propose  an  implicit  solution  method  to  overeome  the  problem. 


MAXWELL’S  VELOCITY  SLIP  BOUNDARY  CONDITION 


In  1879,  James  Clerk  Maxwell  published  a  paper  on  the  viscous  stresses  arising  in  rarefied  gases  [1].  At  the  time, 
a  reviewer  commented  that  it  also  might  be  useful  if  Maxwell  could  use  his  theoretical  findings  to  derive  a  velocity 
boundary  condition  for  rarefied  gas  flows  at  solid  surfaces.  Consequently,  in  an  appendix  to  the  paper,  Maxwell 
proposed  his  now-famous  velocity  slip  boundary  condition.  This  boundary  condition  was  successful  in  predicting 
two  prior  experimental  observations;  (a)  that  a  rarefied  gas  could  slide  over  a  surface,  and  (b)  that  inequalities  in 
temperature  could  give  rise  to  a  force  tending  to  make  the  gas  slide  over  a  surface  from  colder  to  hotter  regions 
(which  had  been  discovered  by  Reynolds,  and  was  known  as  “thermal  transpiration”  —  now  more  commonly  known 
as  “thermal  creep”).  What  has  subsequently  been  overlooked  by  many  current  researchers  is  the  general  form  of  the 
slip  expression  derived  by  Maxwell,  and  this  has  some  substantial  consequences  for  modem  simulations  of,  e.g., 
hypersonic  aerodynamics  and  gas  flows  in  microsystems. 

Maxwell  related  the  tangential  gas  velocity  slip,  w^to  the  tangential  shear  stress,  f ,  and  heat  flux,  q.  Writing 

his  expression  in  vector  form  so  that  it  can  be  easily  applied  to  flows  over  three-dimensional  surfaces,  we  have; 


slip 


up  ^  TP 


(1) 


where  L  =  g  =  Q-(1-  IJ„ ) ,  an  arrow  denotes  a  vector  quantity,  a  is  the  momentum  accommodation 

coefficient  (equal  to  one  for  surfaces  that  reflect  all  incident  molecules  diffusely,  and  zero  for  purely  specular 
reflection),  p  is  the  gas  viscosity  at  the  wall,  A  is  the  molecular  mean  free  path  at  the  wall,  Pr  is  the  Prandtl  number,  / 
is  the  specific  heat  ratio,  p  is  the  gas  pressure  at  the  wall,  is  a  unit  vector  normal  and  away  from  the  wall,  IT  is  the 
stress  tensor  at  the  wall,  1  is  the  identity  tensor  and  2  is  the  heat  flux  vector  at  the  wall.  Maxwell  proposed  a  value 
for  the  slip  coefficient,  Ai=l,  although  a  recent  analytical  solution  of  the  linearized  Boltzmann  equation  at  a  planar 
surface  [2]  estimates  the  “micro-slip”  (the  actual  velocity  slip  at  the  wall)  requires  ^i=V(2/7i)»0. 8.  In  [2],  a  value  of 
^1=1.146  is  also  proposed  when  the  Navier-Stokes  equations  are  used;  this  represents  an  additional  “fictitious” 
velocity  slip  to  accommodate  for  the  Knudsen-layer  stmcture  that  is  not  captured  by  linear  constitutive  relations  [3]. 
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If  the  Navier-Stokes  constitutive  relations  for  stress  and  heat  flux  are  substituted  into  Equation  (1),  expressions 
for  velocity  slip  in  terms  of  flow  gradients  can  be  obtained.  In  his  original  paper,  Maxwell  used  a  one-dimensional 
expression  for  the  shear  stress  (appropriate  for  the  typical  case  he  was  interested  in)  which  made  his  final  result 
generally  applicable  only  to  non-rotating  planar  walls  (i.e.  where  the  streamwise  variation  in  wall-normal  velocity  is 
negligible).  In  scalar  form,  Maxwell  gave; 


^  2>  /Li  dT  (2) 

(7  dn  A  pT  dx 

where  n  is  the  coordinate  normal  to  the  wall,  x  is  the  coordinate  tangential  to  the  wall,  is  the  x-component  of  the 
gas  velocity,  is  the  x-component  of  slip  velocity,  and  p,  T  are  the  density  and  temperature  of  the  gas  at  the  wall, 
respectively. 

It  is  because  of  its  relative  simplicity,  compared  to  Equation  (1),  that  Equation  (2)  is  remembered  as  Maxwell’s 
main  theoretical  result.  However,  for  most  surface  geometries  of  practical  interest,  having  curvature  and/or  rotational 
motion,  it  is  inapplicable  because  it  neglects  that  the  velocity  normal  to  the  wall  can  vary  in  the  streamwise 
direction.  Therefore  a  more  complete  expression  for  the  tangential  shear  stress  is  required  in  Equation  (1).  If 
calculations  straightforwardly  applying  Equation  (2)  are  performed,  they  are  likely  to  miss  important  features  of  the 
rarefied  flow  behavior. 


GEOMETRIC  DEPENDENCE  IN  MAXWELL’S  GENERAL  FORMULATION 

If  the  full  Navier-Stokes  description  of  the  stress  tensor  is  adopted,  for  a  wall  in  two  dimensions.  Equation  (1) 
reduces  to 


y  -  (2-0-)  /aw,  I  'i  p  dT  ^3^ 

a  V  Sk  dx  j  a  pT  dx 

where  Mj  is  the  slip  component  tangential  to  the  wall,  and  u„,  are  the  gas  velocities  normal  and  tangential  to  the 
wall,  respectively. 

The  additional  term  that  features  in  Equation  (3)  but  not  in  Equation  (2)  can  have  a  significant  influence  on  the 
overall  velocity  slip.  For  example,  when  there  is  rotational  wall  motion  (i.e.  wall  motion  in  a  direction  normal  to  the 
surface  with  a  velocity  that  varies  in  the  tangential  direction)  then  even  for  flat  surfaces  (e.g.  a  deflecting  flap)  there 
will  be  a  finite  tangential  velocity  slip.  For  stationary  walls,  surface  curvature  will  also  give  rise  to  a  contribution 
from  the  additional  term. 

Although  the  misapplication  of  Equation  (2)  to  general  geometries  is  widespread,  it  is  not  universal,  and  there  are 
instances  where  curved  boundaries  have  been  treated  appropriately  [4,  5].  Einzel,  Panzer  and  Liu  [6]  derived  a 
boundary  condition  similar  to  Equation  (3)  for  surfaces  with  curvature.  However,  their  boundary  condition  was 
formulated  such  that  slip  due  to  surface  normal  motion  could  not  be  accommodated.  Also,  it  did  not  include  the 
contribution  of  thermal  creep  to  velocity  slip,  nor  was  the  relationship  to  Maxwell’s  general  equation  realized.  Two 
examples  will  show  the  importance  of  maintaining  the  geometric  dependence  that  the  general  form  of  Maxwell’s 
formulation  incorporates. 


Cylindrical  Couette  Flow 

Recent  analytical  and  molecular  dynamics  studies  [5-7]  suggest  that  the  velocity  profile  in  a  rarefied  cylindrical 
Couette  flow  can  become  inverted.  In  the  case  of  a  stationary  outer  cylinder  and  rotating  inner  cylinder,  ‘inverted’ 
means  that  the  radial  velocity  of  the  gas  becomes  greater  further  away  from  the  moving  centre. 

We  have  performed  a  simple  isothermal  calculation  using  a  finite  difference  discretization  of  the  Navier-Stokes 
equations  to  examine  the  influence  of  various  boundary  conditions  on  the  velocity  profile.  The  inner  and  outer 
cylinders  have  radii  of  3%.  and  5?^,  respectively,  and  the  former  has  a  tangential  velocity  approximately  a  third  of  the 
speed  of  sound.  The  gas  is  argon  at  STP  conditions  and  the  accommodation  coefficient,  cj,  is  0.1.  Figure  1  shows  a 
comparison  of  the  velocity  profiles  (non-dimensionalized  by  the  tangential  velocity  of  the  inner  cylinder)  predicted 
using  the  standard  no  slip  condition,  the  conventional  slip  condition  (Equation  2),  Maxwell’s  general  slip  condition 
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(Equation  1),  DSMC  (direct  simulation  Monte  Carlo)  molecular  dynamics  [7],  and  the  analytical  method  of  Einzel, 
Panzer  and  Liu  [6]. 

The  DSMC  method  (being  a  statistical  molecular  dynamics  simulation)  is  often  used  as  an  independent 
numerical  test  in  the  absence  of  experimental  data  [8].  That  DSMC  predicts  an  inverted  velocity  profile  is  strong 
corroborative  evidence  that  the  phenomenon  is  real.  The  conventional  slip  condition,  Equation  (2),  evidently  cannot 
predict  this  behavior;  however,  Maxwell’s  general  slip  condition  produces  just  such  a  velocity.  Close  quantitative 
agreement  between  DSMC  and  simulation  is  not  expected  here,  as  the  degree  of  gas  rarefaction  in  this  problem 
means  that  continuum  fluid  models  are  at  the  limit  of  their  applicability. 


FIGURE  1.  Nondimensional  velocity  profiles  in  simulations  of  cylindrical  Couette  flow.  Comparison  of  no  slip  (■•■), 
conventional  slip  ( — )  and  Maxwell’s  general  slip  ( — )  solutions,  an  analytical  solution  [6]  ( - )  and  DSMC  data  [7]  (O). 


Drag  on  an  unconiined  sphere 

Isothermal  slip  flow  past  an  unconfined  sphere  at  very  low  Reynolds  and  Mach  numbers  was  first  analyzed  by 
Basset  [4]  using  Stokes’  creeping  flow  approximation.  The  effect  of  slip  was  incorporated  into  the  analysis  using  a 
velocity  boundary  condition  for  isothermal  flows  of  an  equivalent  form  to  Maxwell’s  general  boundary  condition. 
Equation  (1).  Basset’s  analysis  showed  that  the  skin  friction  drag,  D„  on  an  unconfined  sphere  of  radius  a,  in  a  flow 
stream  of  velocity  U,  can  be  written  as; 


=  AnjJJa^ 


a(a  -3T)  +  6A 


(4) 


However,  if  the  conventional  form  of  Maxwell’s  boundary  condition.  Equation  (2),  is  used  in  the  derivation  a 
different  expression  for  the  skin  friction  drag  is  obtained; 


Dj  =  4;rjuUa\ 


a(a  +  /i)-2A  '] 
aia-2A)  +  4A  ) 


(5) 


The  disparity  between  the  two  drag  predictions  is  due  to  the  exclusion  of  curvature  effects  from  the  conventional 
boundary  condition.  In  the  limiting  case  of  a  perfectly  smooth  sphere,  such  that  all  incident  molecules  are  reflected 
specularly  (a  =  0),  there  is  no  means  by  which  the  wall  can  transfer  tangential  momentum  to  or  from  the  gas. 
Therefore,  the  drag  due  to  skin-friction  should  be  zero  and,  indeed,  we  find  that  Basset’s  drag  equation  predicts  no 
skin  friction.  However,  Equation  (5)  predicts  a  finite  value  of  negative  skin  friction  drag  (i.e.  a  thrust).  This  non¬ 
physical  prediction  demonstrates  the  importance  of  employing  the  general,  as  opposed  to  the  conventional,  form  of 
Maxwell’s  boundary  condition  for  curved  surfaces. 
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CONSTITUTIVE  DEPENDENCE  IN  MAXWELL’S  GENERAL  FORMULATION 


While  there  is  no  commitment  made  to  the  form  of  the  stress  tensor  or  heat  flux  vector  in  Maxwell’s  general 
condition,  Equation  (1),  normally  the  Navier-Stokes  constitutive  relations  are  assumed.  This  means  that  in  flow 
situations  where  linear  constitutive  relations  are  no  longer  appropriate,  the  generality  of  the  constitutive  dependence 
in  Maxwell’s  expression  needs  to  be  exploited. 

For  example,  consider  the  DSMC  simulation  of  a  planar  Poiseuille  flow  velocity  profile,  U,  shown  in  Figure  2 
[3].  The  function  '-Ailz  IfJ  (the  dashed  line)  is  proportional  to  the  shear  stress,  and  the  function  ‘^y/ldw/dw’  (the 
dashed  and  dotted  line)  is  proportional  to  the  rate  of  strain.  Towards  the  centre  of  the  channel  these  functions 
converge,  indicating  that  a  Navier-Stokes  model  is  appropriate  for  the  “core”  flow.  Near  the  wall,  however,  the 
Knudsen  layer  represents  a  region  of  high  local  non-equilibrium  in  which  the  proportionality  between  stress  and  rate 
of  strain  breaks  down.  The  conventional  form  of  the  slip  boundary  condition  relates  slip  linearly  to  rate  of  strain  so, 
in  this  case  where  the  DSMC  velocity  profile  is  known,  would  predict  a  velocity  slip  (indicated  by  the  triangle) 
almost  three  times  greater  than  the  actual  slip.  In  contrast,  the  general  form  of  Maxwell’s  condition  provides  an 
accurate  prediction  of  the  slip  at  the  wall  (indicated  by  the  circle). 

If  the  flow  field  is  to  be  modeled  throughout  using  the  Navier-Stokes  equations  the  two  functions  shown  in 
Figure  2  are  necessarily  equal.  However,  when  higher-order  continuum  methods  are  used,  such  as  the  Burnett 
equations  [8],  the  simplified  conditions  of  Equation  (2)  and  Equation  (3)  can  no  longer  be  employed  with 
confidence. 


FIGURE  2.  The  velocity  profile  U  and  derived  quantities  -yf  j2r///and,4;/ld«/d«  in  a  DSMC  simulation  of  planar  Poiseuille  flow 
[3].  The  gas  is  argon  and  the  Knudsen  number,  based  on  half  the  channel  width,  is  approximately  0.2.  The  symbols  A  and  O 
denote  the  slip  predicted  by  Maxwell’s  conventional  slip  conditions,  and  the  slip  predicted  by  the  general  form,  respectively.  Here 
the  slip  coefficient  has  been  taken  as  =  V2/ k  ■ 


Thermal-stress  slip  flow 

Thermal-stress  slip  flow  is  a  rarefaction  phenomenon  that  was  originally  predicted  by  Sone  [9].  Using  an 
asymptotic  analysis  of  the  Boltzmann  equation,  he  showed  that  a  tangential  variation  in  the  wall-normal  temperature 
gradient  could  induce  velocity  slip.  Sone’s  configuration  is  a  gas  (initially  stationary)  between  two  stationary  non¬ 
coaxial  cylinders  of  different  uniform  temperature,  Tj  and  T).  In  the  absence  of  thermal  creep  (i.e.  the  boundary 
temperature  jump  is  not  considered)  no  conventional  boundary  condition  has  the  mechanism  to  predict  a  slip  flow. 
Sone,  however,  calculated  the  slip-flow  field  as  shown  by  the  streamlines  and  directional  arrows  in  Figure  3(a).  The 
outer  cylinder  is  held  at  a  higher  temperature  and  this  generates  a  steady  clockwise  circulation  in  the  gas  (and  anti¬ 
clockwise  when  Ti>  Tj). 

For  a  similar  cylindrical  configuration,  and  using  a  finite  volume  code  to  solve  the  Burnett  equations  combined 
with  Maxwell’s  general  condition  (1),  we  obtain  the  steady  flow  pattern  shown  in  Figure  3(b).  In  the  absence  of 
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thermal  creep,  no  flow  could  be  sustained  by  either  of  the  slip  conditions  given  in  Equation  (2)  or  Equation  (3).  It  is 
therefore  clear  that  the  phenomenon  predicted  by  Sone  is  captured  in  our  simulation  because  of  the  inclusion  of  the 
general  form  of  Maxwell’s  boundary  condition. 


FIGURE  3.  Streamlines  of  thermal-stress  slip  flow  between  non-coaxial  cylinders  (uniform  temperatures,  T2>  Ti);  (a)  solution  of 
the  Boltzmann  equation  reproduced  from  [10],  (b)  finite  volume  solution  of  the  Burnett  equations  using  the  general  form  of 
Maxwell’s  boundary  condition.  Part  (a)  reprinted,  with  permission,  from  the  Annual  Review  of  Fluid  Mechanics,  Volume  32 
©2000  by  Annual  Reviews  www.annualreviews.org. 


NUMERICAL  IMPLEMENTATION  OF  THE  SLIP  BOUNDARY  CONDITION 

For  a  Navier-Stokes  simulation  of  a  micro  flow  contained  by  a  planar  isothermal  wall,  the  velocity  slip  is  given 
by 


(6) 

<j  &n 

In  this  most  basic  configuration  there  is  a  physical  requirement  for  the  gradient  of  velocity  to  be  of  the  same  sign  as 
the  slip  velocity.  By  using  a  finite  difference  approximation  for  the  derivative,  this  can  be  interpreted  as  a  numerical 
constraint; 


^<1  (7) 

«/ 

where  is  the  streamwise  velocity  at  the  first  internal  grid  point  a  distance  A«  away  from  the  wall  in  the  normal 
direction.  When  solved  explicitly.  Equation  (6)  is  numerically  unstable  if  the  prediction  of  slip  violates  the  condition 
given  in  Equation  (7).  This  instability  has  been  documented  in  [11]  and  places  a  serious  restriction  on  the 
applicability  of  the  boundary  conditions.  To  demonstrate  this  we  consider  the  discretized  form  of  an  explicit  solution 
to  Equation  (6); 


u=/l- 


An 


JSJl  /  s 

■  = - [m,  -  j 

£ 


(8) 
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where  m*  is  the  explicitly-predicted  slip  and  £=hn/L  is  a  measure  of  the  spatial  resolution  of  the  simulation  (here  we 

take  Ai=a=l).  Assuming  that  the  no-slip  condition  is  the  solution’s  starting  point,  Equation  (7)  can  be  combined 
with  Equation  (6)  to  give  a  condition  for  numerical  stability; 

s  >  Kn  (9) 

In  a  Couette  flow  simulation  with  a  resolution  £■=0.01  (i.e.  100  grid  points  across  the  channel)  the  maximum 
stable  Knudsen  number  would  be  approximately  Kn=0.01.  For  flows  at  the  upper  end  of  the  slip-flow  regime  (up  to 
Kn=0.1)  no  more  than  10  grid  points  across  the  channel  would  be  stable. 

There  are  a  number  of  ways  to  circumvent  this  problem  (such  as  relaxation  of  the  boundary  conditions  [12]),  but 
the  most  efficient  method  is  to  treat  the  boundary  condition  implicitly.  Equation  (6)  is  rewritten  implicitly,  as 
follows; 


An 


A  M,- 
A  +  An 


(10) 


This  form  of  the  boundary  condition  is  unconditionally  stable.  To  show  this.  Equation  (10)  is  combined  with 
Equation  (7)  to  give  the  following  condition  for  stability; 


£•  >  0 


(11) 


which  is  always  satisfied. 


SUMMARY 

In  this  paper  we  have  demonstrated  the  importance  of  using  the  complete  form  of  Maxwell’s  general  boundary 
condition.  Examples  have  been  given  where  conventional  and  simplified  versions  of  Maxwell’s  slip  condition 
neglect  the  geometry-  and  constitutive-dependent  characteristics  of  the  general  formulation  and,  as  a  result,  fail  to 
capture  some  important  physical  phenomena.  We  have  also  discussed  a  common  numerical  instability  associated 
with  the  Maxwell  slip  condition  and  proposed  a  simple  implicit  solution  technique  to  circumvent  this  problem. 
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